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ANALYTICAL AND EXPERIMENTAL STUDY OF COMPRESSIBLE LAMINAR-FLOW HEAT 
TRANSFER AND PRESSURE DROP OF A GAS IN A UNIFORMLY HEATED TUBE 

by Alden F.. Presler 
Lewis Research Center 

SUMMARY 

This is the report of an analytical and an experimental investigation to determine 
what effect acceleration caused by heat input to a compressible gas has on the laminar 
Nusselt numbers and static pressure drops. Helium gas flow in a circular tube was 
studied with uniform entrance conditions and wall heat flux conditions. 

The analytical portion considers the compressible laminar momentum and energy 
boundary-layer equations to be valid throughout the tube cross section and axially down- 
stream from the entrance. A numerical finite difference marching technique is used to 
solve the equations and is adaptable to a wide range of values of physical parameters 
such as heat flux and entrance Mach number. Axial development of velocity and tem- 
perature profiles, Nusselt numbers, and static pressures were calculated for a 200- fold 
range in heat flux. Local bulk Nusselt numbers are little affected by thermal accelera- 
tion, but pressure drops from flow acceleration by heating show large increases due to 
heat flux. 

The experiment closely approximates the geometrical and physical model of the 
analysis. There is a fivefold variation of the experimental heat flux, and entrance Mach 
numbers vary from 0. 035 to 0. 10. Experimentally determined Nusselt numbers and 
overall static pressure drops confirm analytical results, which were obtained for cor- 
responding heat flux parameters and dimensionless tube lengths. 


INTRODUCTION 

This report is a continuation of two previous numerical studies of the effect of heat 
flux level on the behavior of the laminar flow of a compressible gas in a uniformly 
heated circular tube for both heating and cooling. Reference 1 reported the heat- 


transfer and frictional behavior of both adiabatic laminar flow and one diabatic laminar 
case for a helium gas model. Reference 2 examined in detail by numerical analysis the 
flow behavior of helium gas undergoing cooling. It was found in both these reports that 
certain important deviations from classical heat-transfer and frictional results occurred 
because of compressibility effects at low subsonic Mach numbers. 

The present report examines in more detail analytically and experimentally the 
heating effect on local Nusselt numbers and frictional and pressure parameters. In the 
analytical program a 200-fold variation of heat flux level is imposed on the uniform en- 
trance flow. In particular the development of momentum and frictional pressure com- 
ponents in the flow is detailed for this range of heat fluxes. 

An experiment is described in which heat- transfer measurements are made on the 
laminar flow of helium gas in an electrically heated small-diameter Inconel tube. Local 
Nusselt numbers, inlet-to-exit pressure drop, and pressure drop components due to 
momentum change and friction loss are determined for a range of wall- to bulk- 
temperature ratios of 1.4 to 1. 7. These experimental results are compared with nu- 
merical predictions from the analytical model. 

Reference 3 was one of the earliest analyses of laminar flow heat transfer that con- 
sidered temperature variation of the physical properties. This analysis assumed fully 
developed velocity profiles in the circular tube. Constant property analyses for heat 
transfer to laminar incompressible flow for developing profiles are given in references 
4 and 5. 

Laminar flow heat- transfer experiments have not been numerous. A uniform heat 
flux program with both up and down flow in a vertical tube is reported in reference 6. 

The object of that research was to establish criteria for when free- convection effects 
are important relative to the forced convection. An experiment similar to the one re- 
ported herein is described in reference 7. It contains local Nusselt numbers for both 
helium and nitrogen flows and friction factors based on wall shear stresses, but it has 
no data on flow pressure drops. 

Correlation of laminar friction and heat- transfer data of reference 7 is given in 
reference 8. The experimental friction factor based on wall shear stress is normalized 
to the fully developed analytical value by a product transformation of the ratio of wall to 
bulk temperature. 

A recent numerical analysis using the marching technique is reported in refer- 
ence 9. The condition of uniform heat flux was applied to a laminar flow which had ini- 
tially a fully developed velocity profile. Besides the Nusselt number development, only 
friction factors based on wall shear stress were shown. 
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SYMBOLS 


Variables and dimensional quantities: 
a local sound velocity, VyRt 

C specific heat at constant pressure 

hr 

Cy specific heat at constant volume 
D tube diameter 

* 2 

G mass velocity, W/at^. 

h heat-transfer coefficient, q^/ (t^ - tp) 

k thermal conductivity 

L length of test section 

p absolute static pressure 

Q rate of energy dissipation 

q w wall heat flux 

r radius 

t absolute static temperature 

u axial velocity component 

v radial velocity component 

W flow rate, mass per unit time 

z axial distance from entrance 

q viscosity 

p density 

r shear stress 

r average shear stress over axial distance z 
Dimensionless quantities: 

Cp specific heat ratio Cp/C p e 

fp p friction factor based on pressure drop from inlet, with density and velocity 
evaluated at local bulk temperature, ^-^(p e - Pj/^~^Pb u b 

f T p local friction factor based on wall shear stress, with density and velocity eval- 
uated at local bulk temperature, 2r w /p b u‘; 
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Gz Graetz number, [(Re)(Pr) e J/(z/D) 

K thermal conductivity ratio, k/k g 

M local Mach number, u/a 

Nu^ bulk Nusselt number, hD/k^ 

2 

P pressure drop from inlet, (p 0 - p)/p 0 u e 
Pr Prandtl number, C ju/k 

r 

Q + wall heat flux, q w r w /k g t e 

R radius ratio, r/ r w 

Re Reynolds number, GD/ ju 

T static temperature, (1 - t/t e )/Q + 

U axial velocity component, u/u 

V radial velocity component, 2Re g (v/u e ) 

Z axial distance, (z/D)/Re 

Z + axial distance, (z/D)/Re^ 

y specific heat ratio, Cp/C v 

e total emissivity 

9 temperature ratio, (t - t^)/^ - t^) 

p’ viscosity, p/p e 

P’ density, p/p 0 

t' shear stress at wall, t r /p u 

w W W c 

Subscripts: 

b bulk conditions 

t tube centerline conditions 

d downstream mixing box conditions 

e entry (initial) conditions 

FR friction 

i inside tube wall conditions 

j node point in axial position 

k node point in radial position 


4 



M momentum 

o outside tube wall conditions 

w wall conditions 


NUMERICAL ANALYSIS 


The steady-state compressible momentum and energy boundary- layer equations 
and the compressible continuity equation, in axisymmetric cylindrical coordinates, were 
transformed in reference 2 to the following nondimensional system: 


p'(U 9U + v 9U 

az 


dP + A 


SR/ dZ R 9R 


SR 


( 1 ) 


p-cMuii+v^y ( r 1)M * 


sz 


SR 


Q 


= -L/kr^V 

dZ Pr R SR\ SR/ 


q+ \9R/ 


( 2 ) 


M + i l(p'VR)=0 (3) 

SZ R SR 


The dependent variables are the axial and radial velocities U and V, the static tem- 
perature T, and static pressure P. The independent variables are the radial and 
axial distances R and Z. And the physical properties are density p’, specific heat 
Cp, isentropic coefficient y, viscosity ju', and thermal conductivity K. The system 
parameters are Prandtl number at the entrance Pr , entrance Mach number M , and 
uniform wall heat flux Q . 

The transformations reflect the conditions of uniform heat flux, and uniform ini- 
tial flow at the entrance to the tube, as shown in figure 1. 

The transformed perfect gas law is 


1 - yM^P 

P’ = ~ 

1 - Q + T 


( 4 ) 


The appropriate boundary conditions are 
Initial (Z = 0, 0 < R < 1): 
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U = 1 


Centerline (Z > 0, R = 0): 


Wall (Z > 0, R= 1): 



( 6 ) 


( 7 ) 


Thermal and transport properties appropriate to helium gas (ref. 1) are evaluated 
by the special relations 
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The numerical solution of the system described above by the marching technique 
(see appendix) for given Q + and M g , provides values of U(R, Z), T(R, Z), and P(Z) as 
numerical output. Also included as output are the bulk velocity and temperature, U^(Z) 
and T b (Z), respectively. From these, the heat-transfer parameter 


Nu b = 


2 

V T w " T b> 


(9) 


is immediately calculated, where from equation (8) 


K b = “i - ( l - Q \J' 


68 


( 10 ) 


(ref. 1). 

The two friction parameters of interest in this study are obtained from the numer- 
ical data. The first, based on local wall shear, is 


4f .Re, = 
r,b b 


16t w 

^b U b 


( 11 ) 


where 


t r 

t w w 


T = 
W 


^e u e 


The second, based on the overall pressure change, is 


4f. 


pb Re b “ 


2P 

UuZ' 


( 12 ) 


(13) 


where Z + = (z/D)/Re^ is a dimensionless distance based on bulk Reynolds number. 

The overall static pressure drop is the sum of the drop due to friction and the drop 
due to momentum change, 


Ap A Pf fiction + momentum 


(14) 
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o 

Divide these terms by p u , 

c 0 


Ap _ Ap FR | Ap M ^ 

P e “| P e u e Pe u e 

The left side of equation (15) is the dimensionless pressure P. Then the equation can 
be written as 


P “ P FR + P M 


(16) 


where the dimensionless pressure terms on the right have the form (ref. 2) 


4 T 


FR 


w 


D U 

H e e 


(§) 


(17) 


and 


P M = ^b ” 1 


(18) 


The wall shear stress in equation (17) is an average over the length Z. 

The bulk velocities calculated in the numerical program are used in equation (18) to 
calculate the dimensionless momentum pressure drop occurring between the entrance 
and the local axial position. As the overall static pressure P is also calculated as a 
function of axial position, the dimensionless friction pressure Ppj^ is immediately 
calculated from equation (16) for the required local axial positions. 

The relative distribution of the friction and momentum pressure drop contributions 
is immediately given by the ratio P^/Ppp- 


EXPERIMENTAL INVESTIGATION 

The heat- transfer and pressure drop rig, shown schematically in figure 2, consists 
essentially of an electrically heated tube mounted vertically in a vacuum chamber. 

The vacuum chamber is a vertically mounted steel cylinder, 4 feet (1.2 m) high and 
18 inches (0.46 m) in diameter. The interior was sandblasted to reduce outgassing 
time. Access to the tank is through an opening 12 inches (0. 30 m) wide and about 


8 



4 feet (1.2 m) long. A flange has an O-ring seal which holds the working vacuum to 
around 5 microns pressure. The test section is mounted on the inside of the door, 
which is hinged on the bottom, so that when the door is open the test section is outside 
the vacuum chamber and is therefore easily accessible. 


Flow Measurements 

Helium gas from the six-tank manifold is metered through a parallel bank of four 
rotameters (see fig. 2). Gas pressure at the flowmeters is controlled at 30 psig 
(101 000 N/m ). Primary flow control is by means of valves downstream of the rota- 
meters. Further control is by means of a valve on the exit side of the test section. 

The four rotameters were calibrated with helium gas at 45 psia and 80° F. They 

have maximum flow rates of 0.00030, 0.0027, 0.011, and 0.027 pound per minute 
-8 -7 -7 -fi 

(1.1x10" , 1x10" , 3.7x10 and 1x10" kg/sec). This gives a maximum flow rate with 

— fi 

all four tubes at capacity of 0.041 pound per minute (1. 54x10 kg/sec). 


Test Section 

The test section is constructed of Inconel alloy 600 tubing whose outside diameter is 
0. 125 inch (0. 308 cm) and inside diameter is 0. 089 inch (0. 226 cm). Briefly, the test 
section 19 inches (48. 2 cm) long between the heating electrodes and 23 inches (58. 4 cm) 
long between the static pressure taps. 

The general dimensions of the instrumented test section are given in figure 3. The 
plenum chambers for static pressure measurements and the wall thermocouples, are 
discussed in the next two instrumentation sections. 

Inconel was chosen because of its fine high temperature strength properties and 
because its electrical resistance varies but slightly over a wide range of temperatures 
(ref. 10). Inconel alloy 600 has a liquidus temperature of 2600° F (1692 K) and a solidus 
temperature of 2540° F (1669 K) (ref. 11). The physical properties of Inconel alloy 600 
needed for the present study are found in references 10 to 12, and are reproduced in 
figure 4 for the applicable temperature range. 


Pressure Measurements 


The static pressure drop of the gas flowing through the heated test section was ob- 
tained from the two tap positions at the entrance and exit of the tube, shown in figures 
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2 and 3. Each tap position consists of a cylindrical plenum chamber of 3/4-inch 
(1.91-cm) diameter Inconel alloy 600 tubing 3/4- inch (1.91-cm) long and welded to the 
flow tube, with the chamber fed through four 4- mil holes drilled along crossed diam- 
eters . The inside of holes were carefully deburred after drilling. Details of the plenum 
are shown in figure 3. 

Pressure leads of 1/8-inch (0. 318-cm) Inconel tubing are connected to the plenum 
chambers by threaded fittings, and the pressure leads are brought through the vacuum 
chamber bulkhead with nylon fittings. Outside the vacuum chamber, plastic tubing con- 
nects the two pressure leads to a IJ-tube manometer, 4 feet (1.2 m) long, which can be 
positioned at 0°, 30°, 45°, and 60° from the vertical. The manometer is half filled 
with acetylene tetrabromide colored with a red dye. This fluid has a specific gravity of 
2. 967 (20°/4°), has a formula weight of 345. 70, and is safe to use (ref. 13 formula 
e-265, pg. C-308). It serves well as manometer fluid for uses between water, on one 
hand, and mercury, on the other. 

A second U-tube manometer half filled with mercury is attached in parallel with the 
tetrabromide manometer for cases when pressure drops are too large for the first 
manometer. In this case, shutoff valves on the tetrabromide manometer seal it from 
the system. 


Temperature Measurements 

Wall temperatures on the test section are obtained from nineteen 24-gage chromel- 
alumel thermocouples, which were spot welded at the positions indicated on figure 3. 
These thermocouples are positioned alternately 180° around the tube, and are spaced 
closer together at both the exit and entrance sections where the temperature gradients 
are steeper. The thermocouple wires are encased in ceramic beads for insulation and 
thermal shielding. 

Figure 3 indicates a pair of thermocouples spot welded 3 inches apart on both the 
downstream and upstream pressure tap tubing. These thermocouples, installed to 
measure thermal gradients on the plenum tubine for heat loss calculations, are also en- 
cased in ceramic beads. 

The two copper bus bars connecting the power to the Inconel test section within the 
vacuum enclosure also have a pair of chromel-alumel thermocouples spotted on them 
for the same purpose of measuring the thermal gradients. On the upstream bus, the 
couples are placed 8 inches (20. 3 cm) apart, while on the downstream bus the couples 
are 12 inches (30.48 cm) apart. 

Upstream and downstream mixing box temperatures are measured with Inconel 
sheathed chromel-alumel thermocouple probes inserted at right angles to the gas flow 
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in the mixing cans. A sheathed thermocouple probe is similarly placed in the plenum 
chamber of the rotameters in order to measure the flow temperatures for correcting the 
indicated flow. 

All thermocouple leads are brought through the vacuum tank bulkhead in sealed 
couplings and then to switching boxes outside the tank. Temperatures were recorded on 
flight recorder charts and monitored on a self-balancing potentiometer. The mixing 
box temperature difference is also checked with a precision laboratory potentiometer. 


Voltage and Current Measurements 

The current fed to the test section is measured with a precision milliammeter con- 
nected to a current transformer with a step-down ratio of 240:1. Voltages at the six 
taps were measured with calibrated laboratory voltmeters. 


DATA REDUCTION 

The experimental program provides the following output: 

(1) Test section wall temperatures 

(2) Mixing box probe temperatures 

(3) Temperature gradients on the bus bars and the two plenum leads 

(4) Mass flow rate of helium through the test section 

(5) Electrical current in the test section 

(6) Voltage drops across the test section and the graphite electrodes 

(7) Pressure at the entrance to the test section and total pressure drop across the 

working section 

The first concern in data reduction is that of overall heat balance. For this the 
total energy input rate is the sum of the wattage in the test section and the wattage in the 
graphite electrodes. 


^total input ^ test section + ^electrodes 

These wattages are obtained from the measured current and voltages of the appropriate 
taps. 

Energy losses to the environment are the radiation flux from the heated test section 
and the conduction losses in the bus bars, in the plenum leads, and out the two ends of 
the test section. 
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^loss, total ^rad + ^bus + ^plenum + ^tube end 


( 20 ) 


The heat balance then equates the net heat input with the temperature rise between 
the mixing boxes. 


• • • • 

^net “ ^total input " ^total losses _ ^^p^d ” W ^ 

The losses consist of all heat not used to raise the gas temperature, as indicated in 
equation (20). The preceding equality, of course, is rarely achieved in a heat-transfer 
experiment because of the uncertainties in the laboratory meters. However, this cal- 
culation must be performed as a check to give confidence in the calculation of the local 
net heat flux used in later heat-transfer correlations. 

« • • 

In the heat loss calculations, the conduction terms Q bus , Qpienum’ and ^tube end 
are calculated from 


Qcond = 


At 

AL 


( 22 ) 


where the thermal conductivity k is appropriate to the metal under consideration, A is 
the cross sectional area normal to the direction of heat flux, and At/ AL is the tem- 
perature gradient obtained from the pair of thermocouples spaced L units apart (for 
the plenum and bus bars) or the gradient obtained from the test section wall tempera- 
tures. For the present experiment, conduction along the Inconel tubing (i.e. , the test 
section and plenum leads) is no more than 1 percent of the total losses. This is due to 
the low thermal conductivity of Inconel (fig. 4(a)) and the very small cross section of the 
tubing. 

Heat loss due to radiation from the test section to the inner wall is essentially radi- 
ation between concentric cylinders, infinitely long, of which the area of the outer 
cylinder (vacuum tank) is essentially infinite in value compared with the radiating sur- 
face of the small inner cylinder (test section). For this extreme but appropriate case, 
the net radiation exchange for diffuse radiation is given by (ref. 14, pg. 4, eq. (31-9)) 


^rad ec Ttube ~ * enclosure, 


(23) 


where the flux is based on unit area of the outer surface of the tube, the temperatures 
are absolute, a is the Stefan- Boltzmann constant, and e is the emissivity factor of the 
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heated test section. Since the tube temperature varies along the length from inlet to 
outlet, the total radiation heat loss is the integral 

rad = 2jrr o j£ L q rad < 24 > 

along the tube length. We know the tube temperature, however, only at 19 positions on 
the heated portion (fig. 3). We therefore replace the integral by the summation over n 
incremental areas 


^total rad ^^rad,j^ * * * » 


(25) 


where 


^rad,j q rad,j 


= q rad,j ' 2jrr o AL jJ 


(26) 


Here 


Vdj ■ <* 


4 4 _ t * ) 

\ j enclosure/ 


(27) 


where the enclosure temperature was at 311 K (560° R), and t. is the absolute tube 

th J 

temperature at the j thermocouple. The incremental area AA.(= 2 ttt AL.) is taken 

3 o ] 

so that the thermocouple is at the midpoint of the small length AL^ . 

The emissivity of a heated tube is a complicated function of temperature, surface 
conditions (smoothness, degree of oxidation) and previous metals history such as an- 
nealing. Data on emissivity of Inconel alloy 600 is collected in reference 12 and is re- 
produced in figure 4(c). The experimental runs reported herein were conducted at wall 
temperatures less than 867 K (1100° F), and thus the solid lower curve in figure 4(c) 
seems appropriate to this calculation. However, an emissivity value of 0.4 is chosen to 
represent the tube for these data at all tube positions. The choice of this value was 
arrived at from trial heat balances starting with emissivity values 0. 2. The higher 
value of 0. 4 gave good balances for all experimental runs. It appears that the aging 
process, which produces the almost constant and larger emissivities shown in fig- 
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ure 4(c), also affected our test section in its many hours of heating at temperatures of 
800 K and lower. This is probably the physical reason underlying the constant (and em- 
pirically determined) emissivity of 0.4. 

The practical equation for radiation loss to be used for the present data reduction 
is 


^ total rad = 2 >' r o <7£ (^ t f ‘ ‘end ^ ^j) = 2 "V £ £‘j 4 AL j ‘ 2 "V £Lt Ll 


(28) 


where r Q = 0. 125 inch (0. 315 cm); L = 19 inches (48. 3 cm); t enc j = 560° R (311 K); 
e = 0. 4; and a = 5. 73x10“ ^ watts per square centimeter per K^. 

The total radiation loss calculated by the preceding technique, along with the con- 
duction losses from the bus bars, are the total losses Q^ ogg used in computing the net 
heat input of equation (21). 

The local net heat flux needed for computing the local Nusselt numbers are obtained 
from 


^et, j ^ inner, j (^input,j) 77 ( r o “ r i) AL j “ q rad,j AL 


outer, j 


(29) 


for the incremental areas AAi nner (= 7rr. AL.) and AA 0U t er ;(= ^ r 0 AL.) containing 

fh ? J J ? J J 

the j thermocouple. Here the first term on the right hand side of equation (20) is the 
electrical energy dissipated in the tube in the volume contained in the length AL^. Thus 



(30) 


is the expression used to directly compute the heat flux need in the Nusselt number. 

The temperature drop across the thin tube wall is calculated by equation (2. 14-d) of 
reference 15. The most severe drop occurs when the outer surface is perfectly insu- 
lated. In this special case, and for the heat input rates of the present work, the maxi- 
mum temperature drop across the tube wall is in the order of 1 K, which is negligible 
for this study. Consequently, the wall temperature used in the heat transfer calcula- 
tions are the indicated recorder temperatures. 

The local bulk Nusselt numbers may now be computed at each of the j thermo- 
couple positions using the net heat flux from equation (30): 
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( 31 ) 



At each thermocouple position j, the local bulk Reynolds number is, from its definition, 



(32) 


In equations (31) and (32) the values of the helium thermal conductivity k b and viscosity 
P b are determined at the local bulk temperature from the transport property data 
of reference 16, which is also plotted in figure 5. 


RESULTS AND DISCUSSION 
Analytical 

The system of equations and boundary conditions (1) to (8) are solved numerically 
by a marching technique (see appendix) for radial variations of U and T, and for axial 
variations in U, T, and P. From the solutions the Nusselt numbers and friction fac- 
tors, equations (9), (11) and (13), are obtained as a function of position along the tube. 
Calculations are carried out for values of initial Mach number of 0. 0005, 0. 004, and 
0. 06, and of the heat flux parameter Q + = 0. 135, 1. 35, 13. 5 and 27.0. Numerical re- 
sults of axial velocity profiles for several axial positions, initial Mach numbers, and 
heat flux parameters, are shown in figure 6. 

Hornbeck (ref. 17) shows the velocity profile development from initial uniform con- 
ditions for constant property fluids to the fully developed parabolic (Poiseuille) profile 
at values of Z + ~ 0. 1. This orderly profile development becomes distorted for flows 
of compressible fluids (variable density) with heat-induced acceleration. In the entrance 
regions of the tube this distortion is observed in large fluid acceleration near the wall, 
resulting at higher heat inputs in velocity "overshoot, " that is, velocities near the wall 
exceeding the centerline velocity. Comparison of figures 6(a) to (d) shows that the ac- 
celeration overshoot distortion changes from an unobservable contribution at the low 
Q + = 0. 135 to a large contribution around Z + ~ 0. 01 for the high Q + = 27. 0. 

The overshoot distortion is smoothed out in the downstream sections of the tube; 
but the effects of near-wall acceleration remain, which prohibit a true fully developed 
flow (parabolic profile) from being obtained (figs. 6(c) to (d)). 
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The resulting bluntness of the downstream velocity profile, caused by heat- induced 
acceleration, is exactly opposite the predicted heat-transfer effect on the velocity pro- 
files as cited in references 3 and 15 (fig. 9. 16). The pioneer work of reference 3 as- 
sumed a fully developed flow; this condition eliminated transverse velocity V and axial 
velocity gradient 9U/3Z from the momentum equation to be solved. Thus the heat 
transfer affected the flow only through the temperature variation of the fluid viscosity 
and thermal conductivity. For gases this implies a profile that is blunter than parabolic 
for cooling, and a profile that is more peaked (u./u, > 2. 0) than parabolic for heating. 

- ^ • I 

The actual profiles are shown in figure 6(e) for the two cases Q = 1. 35 and Q = -1. 35. 
The cooling data are from reference 2. 

There is a Mach number effect which intrudes only in figure 6(c) and which also has 
a decided flattening effect on the velocity profiles as the local bulk Mach numbers ex- 
ceed 0. 5 and the flow proceeds to choking. For nearly the same Z + the profile at local 
Mach number of 0.67 is much blunter than the profile at Mach number 0. 39. Flow 
changes are occurring very rapidly in this region. No further results for this inter- 
resting effect are available because of limitations in the computational programs at the 
higher subsonic Mach numbers. 

Heating rate effect on velocities near the wall is shown in another way in figure 7, 
where the centerline- to bulk-velocity ratio development with axial distance is plotted 
with Q + as parameter. Isothermal analyses for flow developing from an initial slug 
profile (refs. 1, 18, and 19) are plotted for comparison. The limit of u^/u^ of 2.0 
represents fully developed Poiseuille laminar flow. There is but little deviation from 
the isothermal curve for the low Q + = 0. 135 case. For Q + = 1. 35, profile develop- 
ment follows isothermal up to Z + ^ 0. 01 where the velocity overshoot in the wall region 
begins to reduce the quantity of flow in the center region of the tube, resulting in the re- 
duced u^/u b and the longer approach to Poiseuille limit. There is much greater de- 
viation for Q + = 13. 5 and 27.0, and overshoot effect begins much sooner, at around 
Z + = 0. 001. At the highest heating rates, there is a minimum ratio around Z + ^ 0. 03. 
This means that, while the bulk and centerline velocity both increase under heating, 
relatively more flow is carried in the wall overshoot region than in the center portion of 
the tube in the region 0. 001 < Z + < 0. 005. 

The decrease in the u^/u^ ratio for Q + =13.5 and M e = 0.06 in the region of 
Z + > 0. 10 is due, as mentioned earlier, to the flattening of the profile in the upper sub- 
sonic Mach number regime. 

The inclusion of variable density and inertia- convective terms in the numerical 
boundary-layer heat-transfer analysis has yielded the unusual acceleration effects on 
axial velocity development which are shown in figures 6 and 7. Temperature profiles 
suffer a similar distortion from the same causes. 
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The intensity of acceleration effects on both T(Z,R) and U(Z,R) is directly related 
to the heat flux level Q + , and this result is apparent for the dimensionless tempera- 
tures in figure 8 where there is a progressive 200-fold increase in the wall heat flux. 

The following qualitative discussion of the acceleration effect on the variable 6 is ap- 
plicable to gas flows at all heat flux levels, but it is most graphically illustrated at the 
highest flux in figure 8(d). 

In the entrance portion of the tube, the accelerating gas in the wall region convects 
most of the heat flux from the wall; the central bulk of the flow remains at constant 
temperature up to axial lengths of about Gz~* ~ 0.01. 

Downstream of this Gz ^ position, heat flux penetrates to the centerline of the 
tube. However, the continuing heat-induced acceleration of the gas flow finally produces 
temperature profiles in the usually fully developed tube regions that are blunter than the 
classical results of reference 3. The profiles in the latter work were calculated using 
realistic gas transport properties, but without the disturbing effect of axial velocity 
gradients which are associated with the radial convection that mixes the central bulk of 
the flow. 

The result from reference 3 which is plotted on figure 8(d) for purpose of compari- 
son is at the value of t / 1, =1.7; this corresponds to Gz~* = 0. 07 in the present case 
of Q + = 27.0. 

The deviation of the present results from reference 3 is rather opposite to that 
classical prediction in that the downstream blunt profiles in figure 8 are similar to the 
cooling heat transfer cases of the earlier work. There, with the absence of radial con- 
vection, cold gas at the wall sustained a large radial temperature gradient (due to low 
thermal conductivity), but the interior flow showed a rather blunt temperature profile, 
that is, with shallow curvature, because of the higher thermal conductivity of the warmer 
gas. 

In figure 8(c), for Q + = 13. 5, the additional effect of high subsonic Mach number 
becomes apparent in the difference in the curves for the same Gz~^ = 0. 37 but with 
different initial Mach numbers. The solid curve for Gz~* = 0.37 has a local bulk Mach 
number of about 0.67, the same as the velocity curve in figure 6(d). 

Since wall-to-bulk temperature ratios have been used as parameters in correlating 
high temperature heat transfer data (ref. 8) it is interesting to observe the variation of 
that ratio with axial distance. Figure 9 details the development in ascending values of 
the heat parameter Q + . In all cases, the wall temperatures are close to the bulk tem- 
peratures in the immediate region of the tube entrance, which is due to mixing from the 
very high convection where the wall velocity gradients are very steep, and also to con- 
duction caused by steep local temperature gradients. Farther down the tube, the de- 
veloping velocities give rise to substantial wall-to-bulk temperature differences and, 
hence, increasing temperature ratios. With the approach to a fully developed condition 
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the wall-to-bulk temperature difference becomes almost a constant. However, the si- 
multaneous and approximately equal rise in wall and bulk temperature implies that the 
ratio of absolute temperatures approaches a value of unity as a limit. 

All the numerical cases in figure 9 show a maximum in t w /t^, the value of which 

increases with increasing Q + . The position of the maximum on the G~ 4 axis de- 

. z 

creases with increasing values of the heat parameter Q . This inverse relation be- 
tween Q + and axial position of maximum temperature ratio simply reflects the influence 
of the more rapid rise of t and t^, and the larger t w - t^ difference at higher wall 
heat fluxes. 

Numerical results for the Nusselt number in developing flow are shown in figure 10. 
These cases are all for low values of initial Mach number, so there is no choking or 
near choking effects in these data. The most interesting result of these calculations is 
the insensitivity of the fully developed Nusselt number to the level of wall heat flux. 

Even for values of Q + = 13. 5 and 27. 0, the Nusselt numbers approach the classical 
uniform heat flux Nu^ = 48/11 as Gz~ 4 — °° (ref. 15, p. 239). This result contrasts 
sharply with the cooling heat-transfer behavior of helium in reference 2, where the fully 
developed Nusselt number was about 75 percent higher than classical uniform heat flux 
value. 

The laminar friction parameters of equations (11) and (13) are shown in figures 11 
and 12, respectively, as development with axial distance. The shear friction parameter 
in figure 11 is strongly influenced in the entrance region by the level of the heat flux 
parameter Q + . This is a result of the fluid velocity overshoot near the wall, as in in- 
dicated in figure 7. The steep velocity gradients in the entrance wall region naturally 
give rise to large wall shear. 

The low Q + cases and the adiabatic results from reference 1 approach the fully de- 
veloped Poiseuille value of 4f ^Re^ = 64. The higher Q + cases could approach a 
slightly higher asymptotic value of, perhaps 68. Because of excessive computer time 
needed, these calculations were not carried out quite far enough to definitely establish 
such a value. 

Both the Q + = 13. 5 and 27. 0 cases give identical friction parameter values for 
Z + > 0. 10 when the initial Mach numbers are very low. For a M = 0. 06, however, 
the Q = 13. 5 case exhibits a choking effect when Z > 0. 2. This result was first re- 
ported in reference 1 and is due to the steepening of the velocity profile near the wall 
for the high subsonic flows. This profile development is shown in figure 6(d). 

The friction parameter 4fp ^Re^, based on the pressure drop from the tube en- 
trance to the local Z + position, is shown in figure 12 for the four Q + values. The low 
Q + cases approach the fully developed Poiseuille limit of 4f b Re t> = 64, The hi S her Q + 
cases seem to approach an asymptote of 40, but its more precise value is also not es- 
tablished because the calculations were not carried out for larger Z + because of the 
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excessive computer time required. The pressure drop friction parameter for Q + = 13. 5 
and M g = 0.06 dips sharply for values of Z + > 0.2 which, again, is the region of high 
subsonic Mach number effect. 

The dimensionless pressure drop P used to calculate 4f^ in the previous 

figure, is shown in figure 13 as a development with axial distance Z + and with Q + as 
parameter. For the higher heating case Q + = 13. 5, the friction pressure drop P FR is 
also detailed. This was calculated from equations (16) and (18) in the manner described 
in the theory section of this report. Two cases for Q + = 13. 5 are shown, with 
M g = 0. 0005 and 0. 06. For M g = 0. 06, the choking effect is again seen at the upper 
Z + limits of the calculation. Figure 13 is useful for obtaining the theoretical P values 
needed for comparison with experimental pressure drop data. 

The relative distribution of momentum and friction pressure drops, shown in fig- 
ure 14 as against Z + , is a sensitive indicator of the changes in the velocity 

profiles along the tube axis. This is most true at the higher heat fluxed Q + = 13. 5 
and 27.0. The very steep velocity gradients at tube entrance, as noted previously, 
produce high wall shear stress which overwhelms in magnitude the P M from the small 
rise in bulk velocity in that region. 


Experimental 

Experimental and analytical results for heat transfer and shear friction are com- 
pared in figures 15 and 16, respectively. In figure 15, experimental local Nusselt 
numbers from the present study, and also helium data from reference 7, are plotted 
along with the numerical adiabatic line for developing uniform heat flux case. The first 
10 percent of the tube length are more sensitive to thermocouple errors because the 
wall- to bulk-temperature differences are small, especially relative to the fully de- 
veloped values. Scatter in experimental data is much less pronounced in the downstream 
tube sections, which is demonstrated by the maximum ±25 percent deviation of the data 
around the asymptotic value Nu^ = 48/11 for large Z + (ref. 1). Such deviation is ac- 
ceptable in laminar heat-transfer experiments (ref. 7), and we therefore conclude that 
measurement and data reduction of variables such as mass flow, temperatures, and net 
heat flux, are correct. 

No local wall shear values were obtained in the present experiment. However, 
Davenport (ref. 7) reported local friction factors for his heat-transfer experiments with 
helium gas. His values, based on local bulk conditions, are plotted on figure 16, and 
they are compared with analytical results for adiabatic case (Q + = 0) and for one high 
heat flux case (Q + = 13.5) both of which also appear in figure 11. The experimental 
friction parameters have the most scatter in the entrance region (Z + < 0. 04 for most 
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runs). This is where determination of the axial pressure gradient (needed for calcula- 
tion of the local wall shear stress from pressure tap data) is most inaccurate. In the 
downstream regions of the tube, the experimental friction parameters seem to converge 
in the limit to the analytical prediction. 

Table I gives the pressure drop results for eight experimental runs. Comparison 
is made in the last two columns with analytical predictions corresponding to the same 
Z + value as the experimental. Since the analytical results in figure 13 for Prp and in 
figure 14 for are at values of Q + = 0.135, 1.35, 13.5, and 27.0, these last 


TABLE I. - COMPARISON OF EXPERIMENTAL AND NUMERICAL PRESSURE DROP BEHAVIOR 





Experimental 




Numerical 

Entrance 

Exit 

Entrance 

Wall heat 

Pressure 

Momentum- 

Axial 

Pressure 

Pressure 

Momentum- 

Mach 
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Reynolds number, 

flux, 

drop 

to friction- 

distance, 
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drop 

to friction- 

number, 

number, 

Re e 

Q + 

from 

pressure 

N + 

parameter, 

from 

pressure 

M e 

M o 



inlet, 

ratio, 


4f P ,b Re b 

inlet, 

ratio , 




P 

P M^ P FR 



P 

P m/ P FR 

0.0348 

0.0434 

1890 

0.35 

9.97 

0.061 

0. 186 

68.0 

9.0 

0.080 

.0441 

.0792 

2450 

1.65 

21.8 

. 102 

. 180 

61.9 

20.0 

. 130 

.0459 

.0742 

2440 

1.00 

14.1 

. 107 

.196 

60.6 

13. 5 

. 135 

.0460 

.0781 

2440 

1.29 

19.6 

. 100 

.216 

64. 7 

19.0 

. 130 

.0471 

.0774 

2410 

1.13 

16.4 

. 105 

.206 

62.4 

14.5 

. 132 

.0557 

.0808 

2850 

1.02 

11. 35 

.097 

. 143 

79.3 

10.0 

. 140 

.0622 

.0878 

2880 

1..15 

9.85 

.098 

.127 

81. 5 

8.5 

. 142 

.110 

.201 

2840 

1.05 

12.2 

. 140 

.150 

65.0 

10. 5 

. 140 


two columns in table I are interpolated from figures 13 and 14 using the experimental 
Q + recorded in the fourth column. Agreement between the experimental and analytical 
values of P is quite good, about 10 percent difference in most cases. 

The momentum-to-friction pressure ratio P^j/Ppj^ is more sensitive to experi- 
mental accuracy than is P. Calculation of P^j from equation (18) requires accuracy 
in both P and T^ at the test section exit Z+ to give a reliable U^. Because Pp^ 
is obtained from P and P^ with equation (16), any error in Pj^ is reflected in Ppj^ 
in the opposite direction and is magnified in the ratio P^/^FR* These facts are ob- 
served in comparing the experimentally determined ratios with the analytical-numerical 
P m/ P FR in c °l umn table I. The agreement is not as good as with the static 

pressures P. 

The experimental and analytical investigations complement each other, with the 
analysis guiding the experiment and the results of the analysis ultimately testing the 
former. Thus, it is concluded from comparison of the results of the two investigations 
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that the numerical marching solution of the compressible boundary-layer equations for 
flow in tubes gives realistic development of heat-transfer and pressure functions for the 
subsonic region. 

The experimental runs recorded in table I are for constant heat flux, to a first ap- 
proximation. The second run in the table (M = 0. 0441, Q + = 1.65) does have a 
16-percent decrease in flux level along the tube because of radiation loss from the 
higher wall temperatures, as calculated by equation (30). This results in a Q + change 
from 1. 80 to 1. 52. The Q + = 1.65 in the table I is the averaged heat balance value. 
This run was chosen, because of its more pronounced heat flux variation, to test the ef- 
fects of a flux variation in the numerical program. These results will be compared with 
those using (1) the average constant flux from the heat balance and (2) the initial flux 
which would hold for the entire tube in the absence of radiation loss. 

Table II contains three data columns, the first reproducing the experimental results 


TABLE n. - RESULTS OF EXPERIMENTAL - 
ANALYTICAL COMPARISON OF PRESSURE 
RESULTS FOR VARIABLE HEAT FLUX 


CASE Q + = 1.80 - 1.52 



Experimental 

Numerical 

Q + 

1.65 

1.80 

1.80 - 1. 52 


(Average) 

(Constant) 

(Variable) 

Re e 

2450 

2450 

2450 

Z + 

0.180 

0. 180 

0. 180 

M e 

0.0441 

0.0441 

0.0441 

M (Z + = 0. 180) 

0.0792 

0. 083 

0.083 

P 

21.8 

16.8 

20.0 

U b 

2.97 

3.18 

3.20 

4f pb Re b 

61.9 

58.0 

64.0 

P m/ P FR 

0.102 

0. 128 

0.119 


of the second row in table I, and the last two columns reporting the results of the fixed 
Q + = 1. 80 and the variable Q + calculations. For these flux levels and flux variation, 
there is very little difference in pressure results from the numerical program. The 
agreement with the experiment for this particular run is about the same as using the 
heat balance average Q + =1.65. At Z + = 0. 180 the analyses for Q + = 1.65, 1.80, 
and 1. 80 to 1. 52 all gave the same value of Nu fe (Nu fe = 4. 34). 
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CONCLUSIONS 


The laminar, compressible, viscous flow in a uniformly heated tube was investi- 
gated both by numerical analysis and experiment. The flow of the gas entered the tube 
with uniform velocity and temperature and at low subsonic Mach numbers. The results 
may be summarized as follows: 


Analytical 

1. For any level of uniform wall heat flux, the local bulk Nusselt numbers exhibit 
only small deviations from the constant property analyses. This is in contrast to the 
results for the same initial conditions but with moderate wall cooling for which the 
nominal fully developed Nusselt numbers are about 75 percent higher than the constant 
property analysis value. 

2. Flow acceleration due to heating causes static pressure drops to increase as the 
wall heat flux level is raised. However, the pressure parameter, which is formed from 
the pressure drop, is much less sensitive to heating levels, and in the nominally fully 
developed tube section is about 2/3 of the Poiseuille limit of 64. 0 at the highest heat flux 
of the present program. 

3. Velocity profiles in the fully developed tube section are somewhat blunter than 
Poiseuille profiles because of continued flow acceleration in the region near the heated 
wall. This result is contrary to the sharp profile found in heat-transfer texts where 
profile variations are due only to viscosity decrease from the wall temperature rise. 

4. The overshoot velocity effect near the heated wall in the tube entrance does pro- 
duce large shear stress in this region and contributes to larger local friction factors and 
pressure drops in the entrance than are obtained in constant property or adiabatic com- 
pressible analyses. 

5. The static pressure development down the tube is split into the sum of the mo- 
mentum and total friction pressure drops. The contribution of momentum change rel- 
ative to pressure drop due to friction is larger for the higher wall heat fluxes. 


Experimental 

1. Experimental bulk Nusselt numbers for a fivefold variation in wall heat flux verify 
the analytical prediction of insensitivity of the Nusselt numbers to the heating level. 

2. Static pressure drops from the experimental tests are compared with analytical 
values for the same dimensionless heat flux, and results from the two methods agree to 
within 10 percent. 
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3. The general conformity of the analytical and experimental pressure relation at 
the heat flux levels reported herein are considered as ample justification of the relia- 
bility of results and method for the numerical marching technique analysis of developing 
laminar compressible tube flows with heat transfer. 

4. An experimental run in which the wall heat flux decreases because of radiation 
losses by 16 percent from inlet to exit is compared with numerical calculations with the 
same heat flux distribution. The numerical analysis produces values of the pressure 
parameters P, momentum- to friction-pressure ratio P M //P FR’ and pressure drop 
parameter 4fp ^Re^ which are close to the experimentally determined values of these 
three parameters. 

Lewis Research Center, 

National Aeronautics and Space Administration, 

Cleveland, Ohio, January 22, 1971, 

127-01. 
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APPENDIX - NUMERICAL PROCEDURE FOR FINITE DIFFERENCE 
SOLUTION OF LAMINAR FLOW BOUNDARY LAYER EQUATIONS 


The nondimensionalized momentum, energy, continuity, and state equations (1) 
to (4) are, respectively, 



+ 1 JL (p’VR) = 0 
3Z R 3R 


(Al) 


(A 2) 
(A3) 


1 - yM 2 P 

P’ = — 

1 - Q + T 


(A4) 


The finite difference analogs of equations (Al) and (A2) are obtained by replacing the 
derivatives in these equations by 


au ^ U j+l,k ~ U j,k 
az az 


(A5) 


au ^ u i,k+i ' U j , k- 1 

3R 2 AR 


(A6) 


a 2 u ^ u j,k+i " 2 u j,k + u j,k- i 


3R 


(AR)" 


(A7) 


The mesh grid used for these differences forms are shown in figure 17. It should be 
noted that equations (A5) and (A7) are written in explicit form. With the coefficients 
evaluated at the (j,k) mesh point, the numerical scheme to be solved is 
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Equations (A8) and (A9) are used to calculate the dependent variables U and T 
stepwise at the mesh point (j + l,k) from the variables and their derivatives along the 
radius at the preceding j th axial position. 

The transverse velocity V\ k is calculated from the continuity equation (A3) as the 
integral relation 



(A 10) 


where the integral can be obtained by any appropriate quadrature method. In this work 
the Simpson integration formula is used. 
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The integrand in equation (A10) is expanded using the dimensionless equation of 
state (A4): 


>1 k v, k = -L R fp. « ♦ L Q + « . yU 2 dP\ 

J,K 3,K Jo dZ 1 _q+ t \ az dZ /j 


dR 


(All) 


The derivatives on the left hand side of equations (A8) and (A9) are now expressed 
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At the mesh position (j + l,k) the expressions (A12) and (A13) make use of back- 
ward differences in U and T on the left hand side. The right hand sides make use of 
previous values of these variables, along with their axial and transvers derivatives, at 
position (j,k). 

After the momentum and energy difference equations are solved for U. * and T. , 

3+1 ]+l 

at axial position (j + 1), the pressure gradient at (j + 1) is obtained from the momentum 
equation at the wall (R = 1) by 


= 1) by 

\--l 

ii’R— \] 

H+i 

|r 9R \ 

‘ 9R /Jr=i 


j+1 


= -4 


aft 2 V 


au 

9R 


(A14) 


|r=i 

j+i 


The new value of dP/dZ at the (j+1) axial position is then inserted into the right 
hand sides of equations (A12) and (A13) for the subsequent calculations of an d 

Tfc j + 2 for k = 0 to N. The value of is then obtained from the backward dif- 

ference approximation 



(A15) 


The central difference approximations for the transverse derivatives used within the 
interior mesh points are not sufficiently accurate at the wall (k = N) because the 
k = N + 1 value of U and T are fictitious values within the solid boundary. Conse- 
quently, a different differencing approximation is used for the wall derivatives which 
uses the boundary values and several interior mesh points. This wall difference scheme 
is derived from a Taylor series expansion using backward differences from the wall. 
Thus, for any axial j position, 



(A 16) 
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The N, N - 1, and N - 2 subscripts in (A16) and (A17) refer to the k index. The j 
index has been suppressed for convenience in the derivation of the special wall deriva- 
tives, but it will reappear in the final expression for the pressure gradient where the 
wall derivatives are needed. 

2 2 

Eliminating the second derivative (9 U/3R between equations (A16) and (A17) 
gives the first derivative as 



U 


N- 


■2 " 4U N- 


2 AR 


+ - (AR) 
3 



+ . . . 


(A18) 


or 


dU\ s U N-2 ~ 4U N-1 
9 R/ n 2 AR 


(A 19) 


2 

which is accurate to the order of (AR) . 

Eliminating the first derivatives (9U/9R) N from (A 16) and (A 17) gives the second 
derivative at the wall as 


or 



U N-2- 2V N-l 


(AR)‘ 


+ ar(^\ + . . . 

WL 


9 2 uV ^ U N-2 ~ 2U N-1 
,9R 2 / n (AR) 2 


(A 20) 


(A 21) 


which is of the order AR. 
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By substituting equations (A19) and (A21) into equation (A14), the difference form of 
the pressure gradient equation becomes 


(-) =“ 4 
H + i 


n. U j +1 , N - 2 2U j+l>N-l / , 
M j+1,N 2 + ,Mj 

(AR) 2 
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(A 22) 
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Flow direction 

z=0 

Figure 1. - Model for numerical analysis of uniform heat flux to laminar compressible flow 
in tube with uniform initial conditions. 



Figure 2. - Flow schematic of heat-transfer test stand. 
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Figure 3. - Test stand instrumentation. (All dimensions are in inches (cm). ) 
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Figure 4. - Thermal properties of Inconel alloy 600. 
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Figure 9. - Variation of wall- to bulk-temperature ratio along tube axis. 
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Figure 10. - Numerical results for Nusselt number development with several heat fluxes. 



41 



Friction parameter, 4f pb Re b = 2P/Z + U b 
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Figure 14. - Development of ratio of momentum pressure to friction pressure along tube length 
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Figure 17. - Finite difference mesh for numerical solution of the boundary- 
layer equations. 
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